home *** CD-ROM | disk | FTP | other *** search
/ Collection of Tools & Utilities / Collection of Tools and Utilities.iso / basic / ubmalm.zip / tlehmer.ub < prev    next >
Text File  |  1990-08-22  |  1KB  |  34 lines
  1.    10   ' Driver program for subroutine Lehmer.
  2.    20   print "What prime? ":input P
  3.    30   print "A quad. res.: ":input Q
  4.    40   gosub *Lehmer(P,Q,&R)
  5.    50   print "A solution is: ";R
  6.    60   print "CHECK: ";(R*R)@P
  7.    70   end
  8.  1070   *Lehmer(P,Q,&R)
  9.  1080   ' Method of D. H. Lehmer to solve x^2 = Q (mod P), assuming that
  10.  1090   ' P is an odd prime and Q is a quadratic residue modulo P.  The
  11.  1100   ' result is returned in R.
  12.  1110   ' 5 June 1990
  13.  1120   dim Bit%(400) ' Holds bits of (P+1)\2, handles 120-digit numbers
  14.  1130   local Te,H,W,Vs,Vb,T,Vm
  15.  1140   local I,J
  16.  1150   Te=P@8:if or{Te=3,Te=7} then R=modpow(Q,(P+1)\4,P):return endif
  17.  1160   if Te=5 then T=modpow(Q,(P-1)\4,P)
  18.  1170   :if T=1 then R=modpow(Q,(P+3)\8,P):return endif
  19.  1180   :R=modpow(4*Q,(P+3)\8,P)
  20.  1190   :if odd(R) then R=(R+P)\2 else R=R\2 endif
  21.  1200   :return endif
  22.  1210   H=1:repeat inc H:until kro(H*H-4*Q,P)=-1
  23.  1220   W=(P+1)\2:J=-1
  24.  1230   repeat inc J:W=W\2:Bit%(J)=res until W=0
  25.  1240   Vs=H:Vb=(H*Vs-2*Q)@P:T=Q
  26.  1250   for I=J-1 to 0 step -1
  27.  1260   Te=2*T:Vm=(Vs*Vb-H*T)@P
  28.  1270   Vs=(Vs*Vs-Te)@P:Vb=(Vb*Vb-Q*Te)@P
  29.  1280   T=(T*T)@P
  30.  1290   if Bit%(I) then Vs=Vm:T=T*Q else Vb=Vm endif
  31.  1300   next I
  32.  1310   if odd(Vs) then R=(Vs+P)\2 else R=Vs\2 endif
  33.  1320   return ' End of subroutine Lehmer.
  34.